Separation quality of a geometric ratchet 
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We consider an experimentally relevant model of a geometric ratchet in which particles undergo 
drift and diffusive motion in a two-dimensional periodic array of obstacles, and which is used for 
(N . the continuous separation of particles subject to different forces. The macroscopic drift velocity and 

' diffusion tensor are calculated by a Monte-Carlo simulation and by a master-equation approach, 

D , using the correponding microscopic quantities and the shape of the obstacles as input. We define a 

^"*^ ■ measure of separation quality and investigate its dependence on the applied force and the shape of 

the obstacles. 
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<D 1 I. INTRODUCTION 

During the last decade, ratchets have been the subject of intense research efforts (for a recent comprehensive review 
seeEl). Ratchets are able to produce a directed current of particles although no net average force is applied. Besides 
the fundamental interest in such a somewhat counterintuitive physical phenomenon, their analysis-is important both 
, for the description of natural nonequilibrium transport processes (like "Brownian motors" in cellsotj) and for concrete 
technical applications as rectifiers and separation devices. The various types of ratchets considered so far include 
rocking, flashing and correlation ratchets, where a temporally periodic force, periodic switching of a potential and 
colored non-thermal noise, respectively, induce directed transport in an asymmetric potential au. Apart from these, 
there are "geometric" ratchets, which do not necessarily xeouire any time-dependent forcing but consist, instead, of 
a two-dimensional periodic array of asymmetric obstaclesuTEj, see Fig. Particles are driven by a constant external 
force through the array while they are undergoing diffusive motion. Because of the asymmetry of the obstacles, 

l— ~~ ', the particles' average drift velocity acquires a component perpendicular to the direction of the external force, which 
constitutes the ratchet effect. Since this is. easier to realize experimentally than the time-dependent ratchets, it has 
been proposed and already demonstratedE3 as a device for the separation of charged biomolecules which were subject 
• to an external electric field and underwent diffusive motion in an array of micrometer-sized obstacles produced by a 

£ — . ' lithographic process. 

In this article, we want to take up the latter example and analyze specifically the quality of the separation effect and 
its dependence on various parameters. A similar numerical analysis has been carried out in Ref. [ll], where the ratchet 
effect was investigated for a smooth periodic potential. However, we will emphasize that optimizing the ratchet effect 
alone is not equivalent to optimization of the separation quality. Apart from that, we will discuss the criterion for 
assessment of the separation quality, point out several "trivial" possibilities for optimization and analyze the effective 
change in the diffusion tensor brought about by the presence of the obstacles. The discussion will center around a 
model of particles undergoing drift-diffusive motion on a discrete lattice. This model is analyzed numerically using 
both a Monte-Carlo scheme as well as the numerical solutions of a master equation. 

The remainder of the article is organized as follows: First, we will review briefly the geometric ratchet used for the 
separation of particles, point out the distinction between microscopic and macroscopic drift velocities and diffusion 
tensors and establish a measure of separation quality. Then, we will set up our lattice model. Using the results of 
the numerical simulations, we will present the dependence of the macroscopic diffusion tensor and drift velocity on 
the parameter characterizing the force applied to the particles, as well as on the shape of the obstacles. Finally, we 
will optimize the separation quality for a situation of two particle species that are subject to different forces, for a 
restricted set of obstacle shapes. 
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FIG. 1. The model situation: A particle diffusing and drifting through a periodic array of obstacles. When it has reached 
the final row of the array, its horizontal deflection is registered, which leads to a Gaussian distribution for an ensemble of many 
particles (after coarse- graining over several obstacles). The inset shows the transition rates used in the Monte-Carlo simulations 
on a lattice. 



II. SEPARATION IN A GEOMETRIC RATCHET 

In a "geometric ratchet" , particles drift and diffuse in a periodic potential, where the potential inside each elementary 
cell is asymmetric, leading to a ratchet effect. In the present work, we will only consider the special case of a two- 
dimensional array of impenetrable obstacles. The drift velocity and diffusion constant of the particles are assumed 
to be the same everywhere outside the obstacles. While the diffusion constant is fixed, the drift velocity depends 
not only on the mobility but also on the force applied to the particles. If the latter derives from an electric field, it 
will be proportional to the charge of the particles, which is important for effecting a separation of different particle 
species, in the way it has been demonstrated in Ref. [l^. In this work, differently charged biomolecules have been 
injected at a corner of a periodic array of obstacles. An electric field is applied to the setup so that the particles 
are subject to a force pointing downward. Due to the asymmetry of the obstacles, the average drift velocity has a 
horizontal component. Therefore, after the particles have traversed several rows of obstacles, the center of the particle 
distribution is deflected by a certain amount. When the magnitude of this deflection is sufficiently different for two 
particle species, they may be separated by collecting the particles arriving at a certain row below the injection point. 
This permits a continuous separation of particles, in contrast to electrophoresis. The quality of separation does not 
depend on the strength of the ratchet effect alone but rather on the difference in deflections for two given species. 

A. Macroscopic drift and diffusion 

If the center of the particle distribution moves at a "macroscopic" average drift velocity v, then the slope of the 
line it traces (starting from the point of injection) is given by the ratio v x /v y . The average deflection in the final row 
is obtained by multiplying this ratio with the height of the array. Obviously, if the average drift velocity vector for 
one of the species is just proportional to that of the other one, no separation can result, regardless of the strength of 
the ratchet effect (i.e., the magnitude of the slope) itself. 

Apart from the average drift, the particle distribution will also undergo diffusion, with a "macroscopic diffusion 
tensor" D that will be different from the microscopic one, due to the presence of the obstacles (which may hinder the 
expansion of the particle cloud in some directions, for example) . It is important to know about the diffusive spreading 
of the distribution, since it affects the separation quality. A large difference in deflections for the two species will be 
useless if it is bought at the price of a large width of the respective distributions, which will overlap so that they 
cannot be separated unambiguously. 

The general functional dependence of v and D can be obtained from dimensional analysis. The microscopic param- 
eters entering are the microscopic drift velocity v (pointing along the y-direction), the (isotropic) diffusion constant 



D, the height of the obstacles h, and a collection of parameters describing their shape S (including the aspect ratio). 
The only possibility of forming a dimensionless parameter out of a combination of D, v and h is given by 

(1) 

It is proportional to the microscopic drift velocity v, and therefore to the microscopic mobility fi multiplied by the 
force F. The parameter £ will be used to present the results of our numerical simulations of a lattice model and to 
compare them with the real physical parameters. The most general form of macroscopic drift velocity and diffusion 
tensor is given by: 

V = vv (£, S) 

D = DD (Z,S). (2) 
Here, vo and Dq are dimensionless vector and tensor functions. 

From these considerations, one can already conclude that two species differing only in their mobilities (but not in 
the forces acting on them) will not become separated, in spite of the ratchet effect. This is due to the Einstein relation 
D = /xfcsT, which ensures that the ratio /j,/D will be the same for both species. Since the forces are assumed to be 
the same, v/D is also unchanged, such that the respective values of £ are equal. Therefore, the average drift velocity 
vector only gets scaled when one passes from one species to the other, so that the drift slope remains the same, as 
has been discussed above. 

For small values of the external force (i.e., £), the macroscopic drift velocity depends linearly on the force and, 
therefore, on v (as long as the linear component is not suppressed due to symmetry). Therefore, a "macroscopic 
mobility tensor" p relating v to the force F may be defined, such that v — (p/fx)vn, with ft the direction of the 



microscopic drift. The numerical calculations (see remarks in section III) will confirm that it fulfills the Onsager 
symmetry relations: p is found to be symmetric. Likewise the macroscopic diffusion and mobility tensors are connected 
by the Einstein relation for small values of the external force: D — pksT, where the temperature T is obtained from 
the ratio D/fi of the corresponding microscopic quantities. Physically, this derives from the fact that the equilibrium 
distribution in a setup with a wall at the bottom is given by the Boltzmann distribution, which carries over from the 
microscopic density to the coarse-grained density, whose evolution is governed by D and p. Apart from the fact that 
one has to be still in the regime of linear variation of v for the Einstein relation to make sense, it can only hold as 
long as the force is not so strong as to make the density fall off rapidly over the scale of a single obstacle, because then 
the coarse-graining procedure is no longer justified. This happens approximately at £ ~ 1, which is also the condition 
that has to be reached to see an appreciable separation effect. 



B. Analytical estimate 



The magnitude of the ratchet effect can be estimated analyticallyOiffl for high enough external force (large £). 
Then, one can treat the motion in the direction of the external force as deterministic (neglecting diffusion), so that 
diffusive spreading takes place only in the perpendicular direction. In this simplified picture, the geometric ratchet 
becomes analogous to a time-dependent one-dimensional "flashing" ratchet, with the y-coordinate playing the role 
of time. For a typical obstacle (as it has been used, e.g., in the experiment of van Oudenaarden and Boxenld), the 
diffusing particle distribution (in the shape of half a Gaussian curve) will be split in two parts by the "top part" of the 
obstacle, see Fig. ^. Here and in the following, we assume that there is one connected obstacle. The left part proceeds 
downward further on, while the right part moves one cell to the right. Therefore, the slope v x /v y is given essentially 
by the percentage of particles that have moved to the right in such an "elementary step" , i.e., by an integral over the 
respective part of the Gaussian distribution. For large £ it can be approximated by an exponential (which becomes a 
good approximation if the magnitude of the following exponent exceeds 2): 

v<c\ v x 2b fW 1 / w 2 A 



nw V h \ Ahh 



Here, the dependence on £ has been made explicit, while all the other dimensionless parameters are ratios of lengths 
determining the shape S of the obstacle (see Fig. ||). Obviously, if the force becomes very large, the particles will 
only move down inside "channels", since they do not have time to spread to the left or right. Then the slope oc v x /v v 
becomes very small, as is expressed by this formula. On the other hand, for very small forces neglecting the possibility 
of diffusing backwards in y-direction (or more than one cell in ^-direction) renders this estimate invalid. Qualitatively, 
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however, it is correct that the slope tends to a constant in the limit of vanishing force £ — > 0. Therefore, separation 
is ineffective both at very small forces (since the slopes of two species are the same) and at very large forces (since 
the slopes differ but are both very small). Since the prefactor multiplying £ in the exponent is of the order of one, 
the optimum separation quality for the obstacle shape discussed here will also be reached when £1 and £2 are of the 
order of one. This is confirmed by the numerical analysis below. 




C. Quality of separation 

In the long-time limit, the ensemble of particles of a given species (having started at the injection point) assumes the 
form of a two-dimensional Gaussian distribution, which drifts at a velocity v. Therefore, the distribution of particles 
along the x-coordinate of the final row is also a Gaussian, which is centered around some value (x) and has a width 
a. We will assume that separation is performed by collecting all the particles up to some point xq in one bin and the 
rest in another bin. Ideally, the two bins would only receive particles of a single species (1 or 2). Due to the overlap 
of the two distributions, this is not possible and there is a certain percentage of particles that are attributed to the 
wrong bin. Qualitatively, the optimal choice of xq is one where this "error" is minimized. However, since there are 
two different types of errors (percentage of particles 1 attributed to bin 2 and vice versa) , no unambiguous definition 
of the optimal choice of xq and the corresponding optimal separation quality exists. We suggest to take the separation 
of (xi) and (X2) and divide it by the maximum of the widths <j\ and 02 to arrive at a measure of separation quality 
which is easily evaluated: 

Q = \(xi) - (x 2 )| /max(<ri,£72) • (4) 
For the parameters considered here it seems to be appropriate. 

Another possible definition consists in replacing the maximum by the geometric mean of the two widths: 

Q'=\(x 1 )-(x 2 )\/V^- (5) 

Note, however, that there are situations when this may be a misleading measure, particularly when one of the 
widths a 1 or 02 is much larger than the other. In these cases, Q will probably be better suited. 

Other measures have been used in the literature. For example, in Refs. the authors essentially asked how large 
the relative difference in diffusion constants of two species should be in order to have a separation that exceeds the 
spread of one of the two distributions. 

In any case, any reasonably defined optimal separation quality can only be a function of two dimensionless param- 
eters: the measure Q (or Q') used here and the ratio 01/02. 

The goal is to optimize the separation quality for two given species of particles by varying the applied force and the 
shape of the obstacles (including their size and aspect ratio) . The following parameters are naturally assumed to be 
fixed: the microscopic mobilities and diffusion constants of the species (connected by the Einstein relation for a given 
temperature), the ratio A of the forces (equal to the ratio of charges in an electric field) and the total height H of the 
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periodic array in y-direction. The latter will be dictated by practical considerations (if it could be made arbitrarily 
large, one would do this to get an ideal separation effect). 

We want to argue when and why a nontrivial optimum is to be expected. It has been pointed out above that 
choosing the force too large or too small results in a decrease of separation quality. Likewise, one could make the 
ratchet effect itself arbitrarily strong, by choosing a slanted line of large slope dx/dy (as one big obstacle). However, 
both species would drift along that line and could not be separated either. 

In order to obtain a more quantitative understanding, we use the general functional forms given for v and D, and 
insert them into the measure Q introduced above. More precisely, we use the scaling expressions (x) = Hv x /v y = 
Hxq(^,S) and a = W DH/vao(^, S) for both species, which are assumed to have the same microscopic values of /x 
and D but different forces acting on them, such that t>2 = Xvi = Xv and £2 = A£i = A£. Thus, we obtain: 

Q= \(*i)-(x*)\ 



max(<7i,<72) 
~H VZ\x {£)-x {\ti)\ 



MAO/v%<7 (0) 



max 



(6) 



After the fraction in the last line has been optimized by varying both £ and the shape S of the obstacles, Q could 
be made arbitrarily large by having the height h of an individual obstacle go to zero (at fixed aspect ratio), such 
that the array (of fixed height H) contains more rows. Since £ = vh/ D must remain constant, this means that the 
microscopic drift velocity v has to go to infinity. Physically, the enhancement of separation quality can be understood 
in the following way: Although the slopes will remain the same, the relative size of the diffusive spread decreases like 
\JhjH. Note as well that, due to the same reason, the maximum of the fraction does not exist, strictly speaking. It is 
possible to keep h fixed but still have an effective increase in the number of rows by placing an array of miniaturized 
obstacles inside an "elementary" cell. However, in practice there are obvious restrictions on the force that can be 
applied to the particles, as well as on the minimum size of the obstacles. Therefore, it is only possible to choose the 
force as large as possible and the corresponding value of h in such a way that £ takes on the optimal value under these 
restrictions. 

It is also helpful to consider the behavior of Q for the analytical estimate given in Eq. (||), which is valid for large 
values of £. After H/h rows (and a time t = H/v y ), the average deflection is v x t = v y Pt and the variance of the 
resulting Poisson distribution is 2D xx t = bHP. This gives the relation D xx /v x = 6/2, which is confirmed by the 
numerical results below (for obstacles of type "A" and "B" in Fig. ^). We obtain: 

Q~M (7) 
V b max(/P(e), ^PJXT)) 

By setting P(£) ~ exp(— 7£), with some exponent 7, we can find the approximate behavior of Q: It drops like 
exp(— £7/2) at large £, regardless of A(> 1). On the other hand, the alternative definition Q' given above would only 
decrease like exp(— £7(3 — A)/4). For A > 3, this may even rise at large £, which is probably an indication that, for 
this regime, the expression for Q' is not any longer a good measure of separation quality: The difference in spread of 
the two species grows too fast. 



III. THE LATTICE MODEL 



We have set up a model to study numerically the diffusion of particles under the influence of an external force in 
a periodic array of obstacles: A single particle is positioned on a point of a two-dimensional lattice which consists of 
square fields. During each time-step, the particle either changes its position to one of the neighboring squares (with 
a certain probability) or it remains on its original square. 

The probabilities to move to the right or to the left are both equal to T. In the absence of obstacles, there is 
consequently no net flow of particles in the horizontal direction. The probability to move downward (i.e. in positive 
y-direction) is T + a, the probability to move upward is T — a (see inset of Fig. pi). The probability to remain on the 
original square is therefore equal to 1 — 4r. Because of the different probabilities to move up- and downward, a net 
flow results in the vertical direction. 

The following relations hold for the microscopic drift velocity v in y-direction and the microscopic diffusion constant 

D: 

v = 2a 

D = T. (8) 
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We assume lengths to be measured in lattice constants and time in units of the elementary time-step. For the 
interpretation of the results, only the value of the dimensionless parameter £ is needed: 

£ = vh/D = 2ah/T, (9) 

where h is the height of an obstacle (elementary cell of the array) measured in lattice constants. From this expression, 
it is clear that increasing the spatial resolution (h) means decreasing a/T, such that this ratio vanishes in the physically 
relevant continuum limit. Actually, we have already assumed a <C T in writing down Eq. (g|). Otherwise, one would 
have to take into account that the diffusion resulting from the model stated above is anisotropic, with D xx — T and 
D yy — r — 2a 2 . In our Monte-Carlo simulation, we have usually chosen T to be of the order of 0.1 and a/T to be less 
than about 0.1. 

Obstacles are represented by "forbidden" squares. If the particle would have to move onto such a square, it does 
not move and remains on its original square. All obstacles are arranged periodically. 

Particles start in the top left square of the array. They continue moving until they have reached the final row, where 
their ^-coordinate is saved. Using the results of many runs, the average deflection (x) and the standard deviation a 
of the resulting distribution can be calculated, see Fig. [lj Alternatively, a particle runs for a given number of time 
steps and its final coordinates are registered. In this way, the macroscopic drift velocity v and diffusion tensor D can 
be obtained. 

In order to test the Monte-Carlo simulation (which has been implemented in C++), we have verified that the 
distribution of deflections in the final row is Gaussian (when coarse-grained over several obstacles) and becomes 
independent of the precise starting position after a sufficiently large time (number of rows of the array), and all 
quantities show the correct scaling behavior described by Eq. (g), provided a/T is chosen small enough. The lattice 
resolution (connected with a/T) has been chosen such that the results do not depend on it appreciably any more. 

Furthermore, at small values of the external force (i.e., of v and £), the macroscopic drift velocity depends linearly 
on this force and fulfills the Onsager symmetry relations. Likewise the Einstein relation holds: D = p,ksT. 



IV. MASTER EQUATION 

While every quantity of interest (macroscopic drift velocity and diffusion tensor, average deflection and spread of 
distribution) can be calculated using the Monte-Carlo simulation, it is nevertheless useful to consider a master equation 
solution as well. This is both because reaching a high statistical accuracy requires a large number of Monte-Carlo 
runs and because a discussion of the master equation yields additional physical insights. 

The particle distribution p(x, y) is defined on the lattice, with integer coordinates x and y. In a single time-step, 
the distribution changes by 

Sp(x, y) = T(p(x + l,y)+p(x~l,y)) 

+ (r - a)p(x, y + 1) + (r + a)p(x, y - 1) 

-4Tp(x,y), (10) 

where the quantities on the right-hand side are to be evaluated at time t and Sp = p(t+X)—p(t). The equation shown 
here holds for every site (x, y) which has no neighboring obstacle sites. For each "forbidden" site, the corresponding 
incoming and outgoing rates have to be left out. The temporal continuum limit is obtained by letting T and a tend 
to zero, with their ratio kept fixed. 

At large times, an ensemble of particles that has started at the injection point will be spread over many obstacle 
cells. Viewed on the scale of only a few cells (much less than the total spread), the distribution p is periodic. Therefore, 
we can calculate the average drift velocity by solving for the stationary distribution p defined inside the cell of a single 
obstacle, imposing periodic boundary conditions, v is then given by the average velocity inside the cell, i.e., by 

links 

Here, the sum runs over all links connecting adjacent sites and the "current density" j along each link is obtained 
by multiplying the values of p(x, y) on the two connected sites by the transition probabilities of a jump along the link, 
taking into account the direction of the link and leaving out blocked sites. The distribution p itself is assumed to be 
normalized: ~Y^p{x,y) = 1. 

Actually, every site with no neighboring obstacle sites contributes just p(x, y) v to this sum, where v is the mi- 
croscopic drift velocity. In this sense, the deviation of the macroscopic drift v from v is seen to arise only from the 
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boundaries of the obstacle. This can be understood most easily in the continuum limit, where j = —D\7p + vp {p 
being the continuous distribution). An integration by parts yields: 

v = I jd 2 r = v + D(p pnds. (12) 

J cell Jobst. 

Here, the second integral extends over the boundary of the obstacle and h is the normal vector of that boundary 
(pointing away from the obstacle). Therefore, those boundaries where the density p "piles up" (see Fig. || below) 
contribute most to the change in drift velocity. The second integral can therefore be interpreted as being proportional 
to the "force" due to the obstacle, impeding and deflecting the free flow of particles. 

For the numerical solution, we first set up the (sparse) matrix corresponding to the transition rates on the lattice 
with periodic boundary conditions and then solve for the stationary solution p. This is done by setting p(x, y) = 1 
on an arbitrary non-blocked site and striking out the respective column and row in the homogeneous linear system of 
equations, such that it becomes inhomogeneous and nonsingulaj. We have used a sparse matrix bi-conjugate-gradient 
solver from the LAPACK package of linear algebra routinesEj. For a typical cell size of 128 x 128, the number of 
nonvanishing matrix entries is about 10 5 . In the end, p is normalized and the current density and macroscopic drift 
velocity are calculated. This type of calculation has already been performed by the authors of Ref. |ll|, for a model 
with a smoothly varying periodic potential. 

It is also possible to calculate the macroscopic diffusion tensor D using the master equation. This involves the 
solution of an inhomogeneous linear equation, with the linear operator defined by Eq. (|1C|) and the inhomogeneity 
derived from the solution p of the homogeneous equation. The derivation of the equation and the formula for D can 
be found in the Appendix. To the best of our knowledge, such an analysis has not been carried out before. 

V. NUMERICAL RESULTS 

In this section, we present the results of both the Monte-Carlo simulations and solutions of the master equation. 

In the Monte-Carlo simulation, the relative statistical accuracy of the macroscopic drift velocity (Sv x /v x etc.) is 
given approximately by \J D j (v 2 Nt), while that of D is estimated to be I/a/ZV, where N is the number of samples 
and t the number of time-steps. Independence of the detailed initial conditions is reached when the diffusive spread 
\f~Dt becomes much larger than the width of an obstacle cell (typically 160). We have chosen values of t > 10 7 and 
N > 10 3 in order to fulfill these criteria. 



A. Density distribution and flow field 

The particle distribution p and corresponding current density j resulting from the solution of the master equation 
for a typical obstacle shape are depicted in Fig. [| both for a low and a high value of £. This corresponds to the 
optimal situation for the separation of two species, where one of them is almost not deflected at all (high force, i.e., 
high £), while the other one has an appreciable probability to go to a neighboring cell to the right, due to diffusion 
around the top part of the obstacle. 
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FIG. 3. Distribution p (top), velocity direction field v (middle), and "streamlines" (bottom) for low and high force (left and 
right columns). The strength of the force £ = 12 for the left column and £ = 24 for the right column, with the microscopic drift 
pointing downward. Darker shading signifies increased p. Contour lines for p have been chosen at the same (equidistant) values 
for the two pictures at the top. Note: The streamlines serve to visualize the direction of the current density, while the real 
motion of a single particle is governed by drift and diffusion. The "velocity direction" is the normalized current density vector 
field. Note the stronger pile-up in density at the left "roof" of the obstacle for the stronger force, as well as the distortion of 
streamlines there and in the narrow channel at the right of the pictures. The lattice shown here consists of 160 x 160 sites. 



B. Macroscopic drift velocity and diffusion tensor 



The results for the components of the macroscopic drift velocity are displayed in Fig. |j. These data have been 
obtained for the obstacle of Fig. ||, using both the master equation with periodic boundary conditions and the Monte- 
Carlo simulation. As has been explained in a previous section, the onset of a nonlinear dependence of v y and v x on £ 
at rather low £ implies that the obstacle is well-suited for obtaining a separation effect (rather than merely a ratchet 
effect). 




FIG. 4. Macroscopic drift velocity components v x (circles, left scale) and v y (squares, right scale) in units of D/h, plotted 
vs. £ = vh/D, for the obstacle shape shown in Fig. ^, with a lattice resolution of 160 x 160. The lines indicate the results of the 
master equation, the symbols (which are larger than the error bars) are from Monte-Carlo simulations (t = 2 X 10 7 , TV = 10 3 ). 



In Fig. the components of the (symmetrized) macroscopic diffusion tensor D have been plotted vs. £, for the 
same obstacle. At low values of the driving force, all components of D are generally reduced compared with the 
microscopic diffusion constant D, since the obstacles hinder the free diffusion of particles. At higher forces, the y- 
component is increased, while the x-component D xx is further reduced. This can be understood in the following way: 
The motion proceeds inside vertical "channels" , such that the particle cannot move easily in horizontal (x) direction, 
while the diffusion in y-direction is more or less free. Comparing the data for D xx and v x shows that D xx /v x w 6/2 at 
larger values of £ (to the right of the maximum of v x (£)), as expected: both quantities decrease exponentially. Note 
that the off-diagonal component D xy of the macroscopic diffusion tensor changes sign at about £ ~ 8, which seems 
to approximately coincide with the sign-change in Vt^/Vuj, (see Fig. ^). We have not come up, however, with an 
explanation for the approximate correlation between D xy /D and \/v x fVvy yet. 

The statistical accuracy of the Monte-Carlo results for D is worse than that of the v results, as expected. The 
deviation between the results of the master equation and Monte-Carlo simulation at larger values of £ is reduced when 
the grid resolution is enhanced (i.e., when a/T gets smaller for fixed £). In the examples shown here, a/T takes on a 
maximum value of about 0.1. 









i 1 


i 1 


i 




















o D /D 
A D ID 

xy 

□ D /D 

XX 






















[] 



10 20 30 40 



FIG. 5. Components of macroscopic diffusion tensor vs. £, for the obstacle of Fig. |3[ Circles: D yy /D, squares: D xx /D, 
triangles: D xy /D. Symbols show results of the Monte-Carlo simulation, lines those of the master equation. We attribute the 
deviation at higher values of £ (and, therefore, a/F) to effects of the finite grid resolution. 

C. Optimization of separation quality 

For the obstacle discussed above, the slope v x /v y (proportional to the average deflection (x) in the last row) is shown 
in Fig. [| together with those of other obstacle shapes discussed further below. The exponential decay is consistent 
with the analytical estimate derived in a previous section, see Eq. (||). From the slope of the logarithmic plot, a value 
of about 0.17 has been obtained for the prefactor w 2 /(Ahh 1 ) in the exponent of Eq. (|^), which is roughly consistent 
with the geometrical parameters of the obstacle. Given the slope and the spread a (derived from the components of 
_D), one can obtain the separation quality Q defined in Eq. (||), if one assumes some ratio A = £2/^1 of the forces 
acting on the two species (see Fig. fj]). 
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FIG. 7. Normalized separation quality Q^/h/H vs. £, where £ corresponds to the force on the first species and £2 = A£, 
with A = 2 in this plot. Full lines: results of Monte-Carlo simulation, dotted lines: Master equation. The difference is due to 
the larger values for the diffusion tensor component D xx at high force £ yielded by the master equation. The location of the 
optimum does not change much. See discussion in the text. 



D. Influence of the obstacle shape 



Both the magnitude of the ratchet effect and the separation quality depend very much on the shape of the obstacle. 
We have not performed a systematic search over obstacle shapes for a kind of "global" optimization of the separation 
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quality because of the numerical effort involved. However, there are a few general properties resulting from certain 
geometrical features. These are illustrated by the numerical results for the slope v x /v y and the separation quality Q 
plotted in Figs. || and [?]. They can be summarized as follows: 

The vertical "wall" at one side of the obstacle "A" depicted in Fig. || acts to prevent particles from diffusing back 
to the left, thereby increasing v x and leading to a better ratchet and separation effect. This can be seen by comparing 
against a version with a shorter wall ("B"). The triangular "roof" of the obstacle splits the particle distribution in 
two halves as it drifts downwards. If the external force is high, the particles do not have time to diffuse a sufficient 
distance to the right and will be deflected back by the left side of this "roof" , therefore streaming downward, with no 
net deflection to the right. The horizontal position of the upper tip determines the strength of the force where this 
transition takes place: If it is moved to the left ("C"), a much higher force is necessary. In that case, the slope v x /v y 
falls off more slowly with increasing force £ . At the same time, the value of the slope is generally increased at small 
forces, since more particles are deflected one cell to the right. 

In order to illustrate the difference between having a strong ratchet effect and a good separation effect, we have tried 
a triangular obstacle ("E"), which yields a comparatively large slope that, however, does not change very strongly 
with applied force. By flipping the triaagle along the horizontal axis, an obstacle with reflection symmetry is created 
("D"). This has the peculiar featureOO that, for symmetry reasons, v x is an even function of the microscopic velocity 
v y (driving force £), so the linear mobility at low driving force vanishes. In principle, this nonlinear dependence of 
v x on £ is well-suited for achieving a separation effect. However, it must be kept in mind that at low values of the 
external force the diffusive motion is dominant, so the spread and the overlap of the particle distributions of the two 
species in the final row is significant. Therefore, the separation quality decreases towards £ = 0. 

The separation quality can become zero for a special value of the external force whenever the slope shown in Fig. ^ 
has an extremum as a function of microscopic velocity £, such that two different £ can produce the same slope. This 
occurs for two of the obstacles (D and E in Fig. [?]). 

For most of the obstacle shapes considered here, there is at most a local maximum of the separation quality at low 
forces. The global maximum is expected to occur at much higher values of the force, which may be unattainable in 
the experiment (and are difficult to reach in a numerical simulation with a finite lattice resolution). However, for a 
well-suited obstacle like the one depicted in Fig. [|, the quality peaks at moderate values of the force. 

VI. CONCLUSIONS 

In this work, we have analyzed a geometric ratchet consisting of a two-dimensional array of obstacles, where 
particles perform drift-diffusive motion under the action of a constant external force. We have carried out numerical 
calculations using both a Monte-Carlo simulation and a master equation solution in order to obtain the dependence 
of the "macroscopic" drift velocity and diffusion tensor on external force and obstacle shape. Using these results, we 
have quantified the quality of the separation effect that can be achieved when two species of particles with differing 
microscopic mobilities are injected into the array. Our results show the strong dependence on several features of the 
shape of the obstacles and demonstrate the distinction between a strong ratchet effect and a good separation effect. 
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VIII. APPENDIX: CALCULATING THE MACROSCOPIC DIFFUSION TENSOR USING THE MASTER 

EQUATION 

Solving the master equation for the probability density inside a cell containing a single obstacle allows one to obtain 
easily a numerically exact result for the macroscopic drift velocity v. It is given by the average flow velocity of particles 
inside the cell, i.e., the integral (or sum) over the current density (see Eqs. ( [ill ) and (0)). However, to assess the 
quality of separation of a given geometric ratchet, it is equally important to know the macroscopic diffusion tensor, 
which governs the spreading of the macroscopic particle density (i.e., the density averaged over many obstacles). Its 
evaluation using the Monte-Carlo simulation requires a large number of samples, to obtain a good statistical accuracy. 
Therefore, it is desirable to calculate D using the master equation as well. The steps involved in the derivation of D 
are slightly more involved than the straightforward calculation of v. 

Our strategy is to derive the equation of motion for the macroscopic density p, 
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d t p = D tj didjp - uVp , (13) 

from the analogous equation for the microscopic density p, making use of the slow variation of p. 

If the particle density p is spread over many obstacle cells (as is the case after waiting for a sufficiently long time) , 
it may, in a first approximation, be described by 

P = PoP, (14) 

where po refers to the detailed density that varies on the scale of a single obstacle but is periodic throughout the 
array, po has been obtained before by solving the master equation for a single cell, using periodic boundary conditions 
at the borders of the cell and the restriction for the current density jo to run parallel to the walls of the obstacle. If p 
were constant, this would constitute a (not normalizable) stationary periodic solution to the Fokker-Planck equation. 
However, p is assumed to vary (very slowly), such that this is not a stationary solution and does not fulfill the 
boundary conditions exactly. Therefore, it has to be supplemented by "correction terms", which depend (necessarily 
linearly) on the spatial derivatives of p. Consequently, we take p to be given by: 

P = PoP + 9%P + K^didjp + . . . , (15) 

Here, g % and K % i are as yet unknown periodic functions (that vary on the scale of a single obstacle) . We emphasize 
that, of course, many different microscopic densities p yield the same macroscopic density p. Therefore, p is not 
uniquely specified if p is given. However, in the long-time limit assumed here, p has "equilibrated" and the deviations 
from po are in a one-to-one correspondence with the lowest order spatial derivatives of p. 

Our further strategy is as follows: We will first rederive the known result for the macroscopic drift velocity v, which 
is the constant coefficient appearing in the part of the equation of motion for p that contains the spatial derivatives 
of first order: dtp — — vVp. To this end, we will insert the Ansatz (|l5| ) into the equation of motion for p, keeping 
only the terms containing up to first derivatives of p and eliminating g l using the boundary condition. Using the 
result for v, we will obtain an inhomogeneous linear equation for g 1 itself, which must be solved numerically. The 
relation between g % and the macroscopic diffusion tensor D will be obtained by going to the second order in the spatial 
derivatives of p. Although this involves the unknown function K lJ , we will be able to eliminate it in the same way 
that g l had been eliminated in the first step. 

Let us first derive a boundary condition for g l . The current densities for the first and second term on the right 
hand side of Eq. (jl5|) are given by 

3 = joP - Dp (Vp) 
(ji)j = -DidrfMp) + vig'idip) - Dg%dip, 

(16) 

Here, jo = (-DV + v)po is the current density of po alone. We demand j + j\ to be parallel to the obstacle wall, 
keeping only terms including first order derivatives of p and then canceling these terms. This leads to the following 
boundary condition for g % : 

hj g i = Dhipo , (17) 

where j„t = (— DV + v)g l is the current density related to g l and n is the outer normal vector of the obstacle wall. 

In the next step, we calculate dtp = — V^'+ji + . . .) up to first order in the spatial derivatives of p and demand it to 
equal podtp ~ — pot/Vp, which is essentially the drift term for p. Physically, this equation means that a nonvanishing 
slope of p will lead to an overall increase (or decrease) of the microscopic density p inside an obstacle cell. The detailed 
shape of the distribution within that cell is not changed, only its magnitude. After dropping the overall factor c^p, 
we arrive at 

Vj ff . = -2(J )i + + ^)p • (18) 

It is not necessary to know g 1 in order to obtain v. We integrate both sides of this equation over the whole cell, 
assuming periodic boundary conditions for g l (as well as forpo)- The boundary term resulting from the walls of the 
obstacle contains g\ but it can be transformed using Eq. (|l7|) , such that we end up with an equation where g % has 
been eliminated: 

- D I fiipods = -2 o )id 2 r+ (^ + jjj) . (19) 



13 



The integral on the left-hand side runs along the obstacle wall, while that on the right hand side extends over the 
whole cell. Since 



^ (ja)id f=Vi + Df fiip Q ds , (20) 

we have 

Vi = I {johd 2 r, (21) 



as had been_assumed in the main text already (see Eqs. (|11|) and fllS 

Inserting v into Eq. ( |l8| ) yields an inhomogeneous linear partial differential equation for g % which has to be solved 
numerically (assuming periodicity and the boundary condition Eq. (|l7|)). Note that g l is only determined up to a 
constant multiple of po, since po solves the homogeneous equation. However, as we will see, this does not affect the 
result for the diffusion tensor D to be derived from g l . Further remarks concerning the numerical solution of the 
master equation on the discrete lattice can be found at the end of this appendix. 

The current density related to the part of p that involves second derivatives of p (see Eq. (|l5|)) is equal to 

62)1 = (W) [(-£>V, + Vi)K^ + ..., (22) 

where we have neglected higher derivatives of p. We arrive at a boundary condition for K 13 at the walls of the 
obstacle in the same way as for g l , by demanding j to t = j + ji + J2 to be parallel to the wall. This time, we keep only 
the terms including second derivatives of p. This leads to 



-DV 



K ij = Dn jg l . (23) 



To this order, the time-derivative of p, dtp = — V jtot, includes both the diffusion of p and the drift of the term 
9 l dip~: 

Vjtot = vig'didtp - poDudidip + ... , (24) 

keeping only the second order with respect to the spatial derivatives of p on both sides of the equation. 
As before, we integrate this equation over the cell and use the boundary condition Eq. ( |23| ) at the obstacle walls to 
eliminate K lJ . The resulting expression for D then is given by: 



Dji = D Sji — f hjg t ds + (vj - vj) / g l d V . (25) 



Note that adding Apo to g % (with an arbitrary constant A) does not affect the result for D, due to Eqs. (20) and 
©■ 

For the numerical solution, it is, in principle, possible to discretize the equations ([l7]) and (jl^) for g % as well as 
the expression (25) for D. However, this is guaranteed to coincide with the results of the Monte-Carlo simulation 
only in the continuum limit (where, e.g., ol/T —* 0). In order to have a better agreement even when one is not yet in 
the continuum limit, it is advisable to start directly from the discretized master equation and redo the steps of the 
derivation shown here for the discrete lattice. 

The equation which has actually been solved numerically to arrive at g % is given by: 

Lg> = -p$i - T(p l+ -p*~) - 6 l2 a(p l+ +p t '). (26) 

Here, Lg % corresponds to — Vj'gi . Therefore, L is the matrix kernel which is also used for solving the homogeneous 
equation, Lp — 0, including the same treatment of obstacle walls and periodic boundary conditions (see the right-hand 
side of Eq. ( |l0| ) for the definition of L). p is evaluated at the "current site" (the site which the left hand side refers 
to), while p l ~a.nd p % ~ are evaluated at the neighboring sites, in positive or negative direction i (=1, 2, corresponding 
to x, y), respectively. At obstacle walls, these neighboring sites may turn out to be "forbidden", in which case p l+ or 
p % ~ vanishes. This implements the discrete version of the boundary condition discussed above. (Note that according 
to the convention used here, the microscopic drift velocity is assumed to point in negative y-direction if a is positive). 

In order to evaluate D, we must carry out a sum over all sites at the wall of the obstacle, i.e., those which have 
forbidden sites as neighbors. This sum is denoted by Ylw ^ ne sum extending over all allowed sites in the cell is 
denoted by J2 : 
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\ w J w 



(27) 



The last term vanishes in the continuum limit but is important to ensure that D does not change on adding a 
homogeneous solution Xp to g l . 
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